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O ■ Abstract 
(N 

P3 ' We demonstrate that laser beam collapse in highly nonlinear media can be described, for a large number 

of experimental conditions, by the geometrical optics approximation within high accuracy. Taking into 
account this fact we succeed in constructing analytical solutions of the eikonal equation, which are exact on 
jy^" the beam axis and provide: i) a tirst-principles determination of the self-focusing position, thus replacing 

the widely used empirical Marburger formula, ii) a benchmark solution for numerical simulations, and iii) 

^: 

O , a tool for the experimental determination of the high-order nonlinear susceptibility. Successful comparison 

O ' with several experiments is presented. 

^ ' PACS numbers: 42.65.Jx, 42.15.-i, 42.65.-k 
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Nonlinear light self-focusing is a self-induced modification of the optical properties of a mate- 
rial which leads to beam collapse at a certain point Zgf in the media. This effect, first observed in 
the 1960s, plays nowadays a key role in all scientific and technologycal applications related to the 
propagation of intense light beams yj], like material processing [2], environmental sciences ||3D, 
femtochemistry in solutions [4], macromolecule chromatography [i5|], medicine etc. 



Usually, Zs{ is estimated using the empirical Marburger formula 



atlBH], 



which has been con- 



structed via fitting the results of extensive numerical simulations obtained for the case when the 
refractive index n is a linear function of the electric field intensity n = n{I) = no + n2l, (^2 > 0) 
Is , [l^ 111 , 12 ] . In most modern experiments, however, high beam intensities are used for which the 



IJ. 



linear approximation breaks down, and further contributions to n{I) must be considered 113, l8|. 
For these cases no general mathematical condition for the behavior of and the filament intensity 
has been derived so far. Most theoretical results are based on numerical studies, or on variational 
calculations assuming a fixed beam profile inside the medium (see e.g. Ivl] and Refs. therein). 
An analytical theory, able to accurately describe beam collapse in highly nonlinear optics, is still 
missing. Moreover, it is widely believed that the exact treatment of beam propagation in a highly 
nonlinear medium can only be done numerically ll]. 

In this letter, we construct for the first time analytical solutions for the eikonal equations with 
highly nonlinear forms of the refractive index avoiding any a priori assumptions on the form of the 
beam during propagation. The results obtained are exact on the beam axis within the geometrical 
optics approximation, which we demonstrate to be accurate for many of the situations taking 
place in modem experiments. Our approach permits not only to obtain exact expressions for Zgf 
for different nonlinear functions n(I), but also to find a general mathematical framework which 
corrects traditionally used formulas for the filament intensity calculation [|71]. Since the accuracy of 
the semi-classical approximation can be easily estimated, we can determine and control the error 
in our calculations, which is not possible in the case of the Marburger formula. 

Based on these results we are also able to propose experiments to precisely determine the 
high order nonlinear susceptibility of different materials. Moreover, our results yield a natural 
explanation of the experimentally observed IjjI deviation of the scaling law Zsi ~ /^^/^ for high 
beam intensities. 

We consider the propagation of a linearly polarized laser beam of initially Gaussian shape. 
Starting from the nonlinear wave equation and assuming that the light beam is almost monochro- 
matic and that the envelope varies slowly in space and time, one obtains a generalized nonlinear 
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FIG. 1: (Color online.) Accuracy L = L^x/L^^q of the geometrical optics approximation for different 
media as a function of the laser pulse power Pin- The green, black and blue curves refer to water, air and 
fused silica, respectively. The pulse wavelength is assumed to be 800nm. Pin is given in units of TW for 
air, and of GW for water and fused silica. 



Schrodinger equation (NLSE) of the form 



0. 



(1) 



where £ is the electric field, z the propagation length and ko the wave vector. The second term 
describes wave diffraction on the transverse plane. is the nonlinear refractive index. The 

magnitude of the contributions of the diffraction and the nonlinear effects to the beam propagation 
can be estimated through the comparison of the characteristic distances Ldig and L^i, at which the 
beam suffers considerable changes 171]. Then, L = L^i/ L^is is a messure of the accuracy of the 
geometrical optics approximation: if L ^ 1, diffraction can be neglected. 

The main contribution to is usually given by the Kerr cubic term n2\£\'^. Therefore it 

is natural to define a nonlinear length Lni = 1/(A;o^2-^o)j where Iq is the intensity of the beam at 
the entry plane of the nonlinear medium. The diffraction length is deffined as Ldifr = no/coWg/S, 
where wq is the initial beam radius In Fig {T] we plot L as a function of the initial beam power 
-Pin (-Pin = -^0^0^/2) for different media. In many recent experiments L ~ 0.05 or smaller (see flTO 
and Refs. therein), thus making the geometrical optics approximation valid. 

We represent the electric field £ in Eq. ([T]) in the eikonal form: £ = \^ exp{ikoS). Neglecting 
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the diffraction term we obtain 

d.S = -^{d^^Sf + nil), dj = -dMd.S). (2) 

As a next step we rewrite these equations in a dimensionless form by considering x and z in 
units of wq and the beam intensity / in units of Iq. We also introduce the dimensionless variable 
V = dxS and differentiate the first of Eqs. Q with respect to x. The refractive index term yields 
dxn{I) = aipdxl, where cp = ip{I) = din{I) and a = n2/o is a dimensionless parameter. The 
order of magnitude of a in a large number of modern experiments lies below 10~^ Ol]. It is 
therefore reasonable to consider a as a small parameter. 

The boundary value problem to be solved can be summarized as 

dzV + vdxV — aifdxl = 0, dzl + vdxl + IdxV = 0, 

v{0, x) = 0, J(0, x) = exp{-x^), (3) 

which describes the propagation of an initially coUimated Gaussian beam in a nonlinear medium. 
Solutions of Eqs. ([3]) and their derivatives can exhibit singularities for particular values of z. Ana- 



lyzing these points, we obtain the nonlinear self-focusing position Zsf of the laser beam lIlSll . 



Following Ref. UM, one can notice that the system ^ is linear with respect^ the first or- 
der derivatives. Therefore, it is convenient to use a hodograph transformation [16] in order to 
transform it into a linear system of the form 

dwT - -dix = 0, + (^diT = 0, (4) 

where r = Iz, x = x — vz, w = v/a. The boundary conditions are transformed as: for w = 0, 



X = Vl^(V^) andr = 0. 



li 



Taking into account the smallness of a we solve Eqs. dH) by proceeding in two steps [llTI . 
First, accounting for the boundary conditions, from the first of Eqs. (H)), we find r as a function of 
X, / and w: 

r = -w/{2x^). (5) 

Then, substituting Eq. ([5]) into the second of Eqs. dH), we obtain a closed partial differential equa- 
tion for the variable x 

r. aw aw 

OwX + TT^dix + - — i^Oiif = 0. (6) 



Integration of Eq. ^ results in two invariants 

X^ = '^i, J ^-'dl -ar^ = ^2- (7) 

Then, the final steps for the construction of the desired solutions of Eqs. ([3]) can be summarized 
as follows. With the help of Eqs. (|7]) we express / and x as functions of the integration invariants: 
/ = /(\E'i, \E'2) and x = ^2)- Then, we require that according to the boundary conditions 

for r = the equation I(^i,'^2) = exp[— x(\l'i, \1'2)^] must be fulfilled. Finally, substituting 
\E'i and \E'2 into the relation above and returning to the original variables, we get the solution of 
Eqs. ©. 

The scheme presented above allows us to find analytical solutions of the optics equations for 
different types of nonlinearities ip{I) . For high field intensity, the refractive index of most materials 
contains nonlinear contributions additional to the Kerr-term ^2/. Usually, they are modeled as a 



power function of the intensity in the genera 



form pi . Physically, this term can be attributed to 



the fifth order nonlinear susceptibility n^P 11^13] or to the material ionization cxkI^'^ , where K 
is the number of photons absorbed, and ax being the multiphoton ionization (MPI) cross section 

m. 

In the following we use our approach to determine under which conditions self-focussing takes 
place or can be prevented by the high order nonlinearity, and arrive at a new general mathematical 
equation to obtain the filament intensity, which improves previous theories. Let us first consider 
a system having a nonlinear part of the refractive index of the form n{I) = 722/ — n^P (i.e., 
K = 2). In this case cp = 1 — [31, where (3 = 2niIo/n2. Substituting ip into Eqs. (|7]), we get 
^2 = -l//31n(l - (31) - aP and ^1 = x(l - f3I). Thus, the solution to Eqs. © reads 

V = -2alzx{l - - 2alz^{l - 131)]. 

After differentiating these expressions with respect to x and z, solving the obtained system of four 
algebraic equations with respect to d^I, dtv etc., and substituting the resulting expressions into 
Eqs . Q, one can verify that the obtained solutions are exact on the beam axis (x = 0, f |x=o = 0) 



. The on-axial beam intensity distribution is given by 



aPz'^(3 = In 



(8) 



Analyzing the Eq. ([8]) we find a critical value Pc ~ 0.175. For Pc > 0.175, Eq. ([8]) has no 
special points, the on-axial intensity monotonically increases and reaches a saturation value /gat- 



By analyzing the asymptotic behaviour of J = I{z) we obtain /gat = Note that this value 
fuUfils the condition 1 — /3/sat = V^(-^sat) = 0. For Pc < 0.175 there is an interval [zi, Z2] on the 
beam axis where the solution I(z) is not unique. The first point zi can correspond to development 
of an instability in the beam, while the second point Z2 corresponds to z^f, where the intensity 
sharply increases, and the beam is compressed into a single filament with /gat- 

For materials described by n(I) = 722/ — rieP, we have K = 3, = 1 — pP, (3 = 3n6/o/n2 
and the on-axial intensity distribution is given by 

arctanh(-y/5/) — aPz^\J^ = arctanhy^. (9) 

From the analysis of the asymptotic behavior {z 00) we obtain [3^ ~ 0.05, I^sx = By 
inspection, we realize again that, as for K = 2, the intensity of the beam saturates when = 0. 

ues of /sat for K = 2 and 3 obtained here are different from previous theoret- 



Notice that the va 

ical estimates [l7L[l9L|22D, which were obtained assuming that the intensity in the filament saturates 
when the nonlinear terms in n(I) compensate each other From the present results we see, 
however, that this is not the case. Upon propagation, the beam tends to reach the on-axial value 
of the intensity which maximizes the index of refraction at the beam axis. In other words, not the 
nonlinear refractive index, but its variation should be zero: 

dMi)\j^^, = 0. (10) 

This new condition is general, independent of the medium or material and represents one of the 
predictions of this letter, which should serve as a basis for future calculations. Note that such a 
general mathematical condition would have been impossible to obtain on the basis of numerical 
simulations. 

We now apply our theoretical scheme to study the concrete problem of femtosecond laser pulse 
propagation in air, which is relevant due to a large number of applications and whose description 
is still a subject of discussion (see, e. g. Refs. O, ISQ and Refs. therein). The nonlinear refractive 
index of air is taken in the following widely used form 

n = n2[I + /(/)] + n,[P + f{P)] - (11) 

where the first term describes the Kerr response involving a delayed (Raman) contribution /(/) 



K 



^ I-oo e-^*-^'^l^'<I{t')dt'. na and tk are known to be equal to 3.2 x lO^^^ cm/W^ and 70 fs [12311, 



respectively. 77.4 = x*^^Y2^o> where x*-^-* is the fifth order nonlinear susceptibility its exact value is a 



subject of some controversy Isiiy]. The most accepted estimates lie around ~ 10~^^ cm^/W^ [7]. 
In the last term of Eq. (fTTI) . p(/) refers to the density of free electrons and pc denotes the critical 
density above which the plasma becomes opaque. A rough estimate yields p{I) ~ ctxI^ Puttp, 
where pat is the atom density pat = 2 x 10^^ cm^'^ and tp the pulse duration. K = 8 for the MPI 
with a pulse of 800 nm, and (Xg = 3.7 x lO^^^ cm^VWVs OD. 

A numerical solution of the NLSE for the pulse with Pin ~ 0.08TW and wq = 3mm with 
n{I) given by Eq. ^U}) with ^4 = gives Zs{,i = 128.18cm Ref JM. Note that if n{I) = ^2/, 



the self-focusing distance is given by the Kovalev formula [IllL 11211 Zgf = Wq/ (2^722/0), which 
is exact under the geometrical optics approximation and for initially Gaussian beam shape. For 
the experimental conditions of Ref. [24!] this yields f = 127.66cm. This confirms our initial 
statement that for many experiments the diffraction (and in this case also the plasma defocusing) 
can be neglected. 

Recently reported experimental results on air clearly indicate that Zgf scales as 1/ for a 
relatively low initial pulse power Pin. However, for powers above 500 GW (7n ~ 5 x lO^^W/cm^) 
a qualitative change is observed and Zst depends on the power as ~ 1/ Pin \{l4\ . This behavior was 
attributed by authors to noise effects in the beam. On the other hand one can notice that the value of 
500 GW in the experiment of Ref. u4\ corresponds to the case where the terms 722/0 = 1-7 x 10~^ 
and u^Iq ~ 3. 1 x 10^ in the refractive index ([TT]) become comparable. Therefore, we believe that 
the change in the power dependence is mainly driven by the contribution of the highly nonlinear 
terms (fifth order susceptibility). Moreover, from the analysis of the experiment of Ref. nlm within 
our theoretical scheme we can draw important conclusions regarding the form of the nonlinear 
refractive index. 

n _ n 

Using the value of Jq given in Ref. [114'], and assuming 724 ~ =1=10 cm^/W^ [7], we get 
P ~ ±0.35. Comparing this value with Pc ~ 0.175, we see that, if 77,4 < 0, then no self-focusing 
would be observed. Since self- focusing is indeed observed, we conclude that > 0. The behavior 
of Zsi as a function of the initial beam intensity obtained from our theory is shown in Fig. [2] and is 
qualitatively compared to the measured z^f. Note that one could use our results to find an accurate 
value of n4 by performing a similar experiment with a controlled initial Gaussian beam profile 
(i.e., without noise) and fitting the measured curve Zsf(/o) to the solution of Eq. ([8]). 

Based on the results of this work it is clear that one is able to determine experimentally the 
nonlinear optical constants by measuring /sat- As we discussed above, the last term in Eq. ([TT]) 
does not change Zgf significantly. However, its value is crucial for /sat- In Ref. [9] 2;sf and /sat 
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FIG. 2: (Color online.) Dependence of the self-focusing position on the beam intensity. The blue curve 
refers to the dependence l/\/7 (obtained for low intensities), whereas the black curve shows the deviation 
for high intensities due to influence of the fifth-order nonlinearity = 10~^^cm^/W^ . Inset: Experimental 



points for Zfs(Pin) from Ref. 11411 



for air were measured. In the experimental setup a coUimated beam with FWHM diameter of 
d ~ 4 mm (wq = d/V2 ln2), a pulse duration (FWHM) of 450 fs and a pulse energy Ei^ ~ 20 mJ 
was applied. We took these values and used them to calculate the on-axial intensity distribution 
from our analytic equations. The comparison between our theory and experimental results of 
Ref. [25] is shown in Fig^S] The temporal pulse compression due to MPI is estimated along the 
similar lines with Ref. 1231]. For the refractive index given by Eq. (fTTI) . we calculate /sat on the 
basis of Eq. (fTOl ). The value of chosen in a way to provide a satisfactory fit of experimental 
date for both Zst and /sat is equal to 8 x 10^^^ cm'^/W^. It is necessary to point out that for a more 
accurate determination additional experimental studies are required. 



In Ref. n25\\ an additional experiment was reported using a shorter pulse of duration (FWHM) 
50 fs and energy 19 mJ. For the same model input-parameters, we obtain: /sat = 1-26 W/cm^ 
and the fluence F = 0.67 J/cm^ which are in a good agreement with the experimental results 
/sat = 1-3 W/cm^ and F ~ 0.6 J/cm^, correspondingly. 

In Table U we demonstrate the influence of high order nonlinearities on the self-focusing po- 

,-1/2 



sition. The Marburger Formula (MF) reads z,f = 0.367^o ((^ Pm/i^cr - 0.852)2 _ o.0219j 
here zq = tcwq/X ~ 45m, P^r = \^ /'2,iin2 — 3.3GW. The "effective critical power" mentioned in 



the table means a replacement Per Per' where the relation P^ = 2Pcr was used [|24ll . From the 



Table H] one can see that other analytical approaches which do not account for high order nonlin- 
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FIG. 3: Filamentation in air. On-axial fiuence versus the propagation distance. Circles are the experimental 



results from Ref. 1251]. The black curve represents our analytical solution. 



TABLE I: Comparison of the predictions for in air calculated using different theoretical approaches and 
compared with the experimental data of Ref. [25']. 



Experiment 


~ 8 m 


Marburger formula (MF) [9J 


6.1 m 


MF with an effective critical power [7] 


9.9 m 


Kovalev formula (714 = 0) il 11 


9 m 


Eq. dH) for = 8 X lO-^^cm-^/W^ 


8.3 m 



earities cannot provide a satisfactory prediction of Zgf in modem high-intensity experiments. 

Finally, and for the sake of completeness, we present results obtained by applying our theory 
to other types of nonlinearities in order to predict the behavior of self-focusing for other cases of 
current physical interest il3ll . If the nonlinear refractive index has the form n{I) = a//(l + (31), 
the on-axial intensity distribution is given by the expression 

+ -aSpI^z^]'/^ -1 = (3. 
For n(I) = 1 — e~^^^^, the on-axial intensity is given by 

e^' - az'p'l' = 
and for a polynomial form n{I) = (k 2) by 

I^[r^ + k{k - 2)az^] = 1. 
No intensity saturation has been observed. 

9 



The singularities in solutions of equations above correspond to the Zgf for the given n{I). It 
should be noted that the approach suggested in this letter can be generalized to arbitrary n(I). In 
general, the solutions can be found semi-analytically by interpolation of the integrals in Eqs. ([7]). 

Summarizing, exact solutions for the self-focusing length and the filament intensity /sat for 
several different forms of the nonlinear refractive index in the framework of geometrical optics 
were obtained. Depending on the experimental conditions, these solutions can be very accurate, 
describe the essential physics of the problem and explain different independent measurements. 
The analytical expressions obtained for the dependence / = I{z) constitute a clear improvement 
with respect to the empirical Marburger formula. 

This work has been supported by the DFG through SPP 1 134. 
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